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Abstract 

Starting from the hypothesis of a second order transition we have 
studied modifications of the original Heisenberg antiferromagnet on a 
stacked triangular lattice (STA-model) by the Monte Carlo technique. 
The change is a local constraint restricting the spins at the corners of 
selected triangles to add up to zero without stopping them from mov- 
ing freely (STAR-model). We have studied also the closely related 
dihedral and trihedral models which can be classified as Stiefel mod- 
els. We have found indications of a first order transition for all three 
modified models instead of a universal critical behavior. This is in ac- 
cordance with the renormalization group investigations but disagrees 
with the Monte Carlo simulations of the original STA-model favoring 
a new universality class. For the corresponding x-y antiferromagnet 
studied before, the second order nature of the transition could also 
not be confirmed. 
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1 INTRODUCTION 

The critical properties of frustrated spin systems are still under discussion 
[|]. In particular no consensus exists about the nature of the phase transition 
of an antiferromagnet on a stacked triangular lattice with nearest neighbor 
interactions. The symmetry group for a Heisenberg antiferromagnet at high 
temperatures SO3 changes to SO2 for the usual antiferromagnet, but the 
symmetry is completely broken in the low temperature phase of the frustrated 
antiferromagnet. This difference in the breakdown of symmetry between the 
frustrated and non frustrated cases should lead to different classes of critical 
behavior. The controversial point for the stacked triangular antiferromagnet 
(STA) is whether the phase transition is of second order with a new chiral 
universality class [fjj or whether its true nature is a weak first order change. 

Monte Carlo studies || f|, |^, 0, § favor a new universality class whereas 
renormalization group studies indicate a first order transition, since no stable 
fixed point can be found for the Heisenberg case in order e 2 with e = 4 — d 
(d is the dimension of space) §| and also for an expansion up to three loops 



in d = 3 |TIJ . The results for the critical exponents of the numerous Monte 
Carlo simulations are listed in table |l|. They agree reasonably well and also 
with the experimental values in table [|. One notices that the new critical 
indices are quite distinct from the standard ones for the Heisenberg model 
with symmetry breaking SO3/SO2 listed in table 0. The only flaw is the 
small negative values of 77, since this exponent should always be positive 
0,0- 



The discrepancy between the renormalization group analysis and the 
Monte Carlo simulation can be removed by assuming that the first order 
phase transition occurs with a correlation length larger than the diameter 
of the cluster studied. Since this size cannot be substantially increased with 
the present technical means and skills we have tried to investigate modified 
lattice spin systems which should belong to the same universality class. The 
idea is to shorten the correlation length by these modifications in order to 
reveal the first order nature of the transition. This strategy was successful 



for the x-y STA-model we investigated before |13| and we could show that 



in fact the phase transition is of first order in this case. 
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Following Zumbach one can analyze such a weak first order transition 
in terms of the renormalization group approach as due to a fixed point in 
the complex parameter space (More precisely it is only necessary to have a 
minimum in the RG flow). A basin of attraction for the real parameters will 
be generated and "mimic" a second order transition with slightly changed 
scaling relations. The unphysical negative values of rj in table [I] could be 
corrected using Zumbach's approach. 

Investigations about the phase transitions of the STA-systems and similar 
helimagnetic systems 0, || are not only of interest for the field of magnetism 
including the experimental studies. Phase transitions in superfluids, in type 
II superconductors, and in smectic-A liquid crystals should be similar in 
nature to the ones of frustrated Heisenberg magnets. 

The numerical simulations are closely connected to a similar study for 
the x-y spins [O. We study the classical Heisenberg spin system on the 



stacked triangular lattice by fixing for selected triangles the spin direction to 
a 120°-order which would correspond to the ideal antiferromagnetic order on 
a triangular lattice as in Fig. |3|. The common orientation of a 120°-cluster 
is still free. The principle behind this construction is that modes removed 
by the 120°-rigidity are "irrelevant" close to the critical temperature. What 
is relevant is the common direction and orientation. The antiferromagnet 
STA and the rigid antiferromagnet STAR should therefore be in the same 
universality class. These considerations can also be applied to Stiefel's 



or dihedral model |15| . A third model studied is the right handed trihedral 
model. It is constructed by adding a third vector given by the vector product 
of the two of the dihedron. As we do not add new degrees of freedom the 
model should again belong to the same universality class. 

That the behavior close to the phase transition should not change if one 
modifies the models by constraints can only be expected from a second order 
transition. However, if one assumes that the transition is of first order, but 
not visible since the barriers between the two phases are too low for the spin 
clusters one can simulate, then building in constraints in the models and 
simulations should make a difference since the barriers are more difficult to 
overcome. 

In section II we present a review of the RG studies. In the following 
section the three models studied are presented. The details of the simulations 
and the finite size scaling analysis are described in section IV. Results are 
given in section V. The discussion is in section VI where we will review also 
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the experimental studies and the other RG studies in 2 + e expansions. 



2 Renormalization group studies and com- 
plex fixed points 

2.1 Renormalization group studies in d = 3 and d = 4— e 

The renormalization group studies have motivated for the re-analysis of 
the frustrated Heisenberg antiferromagnet with the Monte Carlo technique. 
Therefore we repeat briefly the main points. For more details we refer to 
Antonenko and Sokolov [TIJ . 



The Hamiltonian one considers is 



H = - I d'\r 



with cf) a a complex vector order parameter field. Summation of repeated 
indices a, (3 = 1, . . . N is implied with N = 3 for the Heisenberg case. The 
deviation from the transition temperature is given by r$ and there are two 
coupling constants «o and Wq, both positive for a non-collinear ground state, 
that is the two vectors forming the complex vector should not be parallel. 
For the non frustrated case there is only one field and one coupling constant, 
but for the superfluid 3 He the above complexity also arises |16|, [I7| . 



Renormalization group calculations for d = 4 — e [16, 18, EJ or also directly 



for d = 3 [10] show that four fixed points exist: 



1. The Gaussian fixed G point at u* = v * = with mean field critical 
exponents. 

2. The 0{2N) fixed point H at w* = and u* = u H ^ with 0{2N) 
exponents (see table 

3. Two fixed points F + and F_ at location (uf + ,wf + ) and (uf_,Wf_) 
different from zero. These are the fixed points associated with a new 
"crural" universality class. 

The existence and stability of the fixed points depend on the number of 
components iV: 
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N > N c : four fixed points are present but three are unstable (G, H, 
F-) and the stable one is F + . Therefore the transition belongs to a new 
universality class different from the standard O(N) class. If the initial 
point for the RG flow is to the left of the line (G, F_), see fig. ||a. the 
flow is unstable and the transition will be of first order. 

N = N c : the fixed points F_ and F + coalesce to a marginally stable 
fixed point. One would think that the transition is "tri critical" but 
the exponents are different and not given by the tricritical mean field 
values. The reason is that there are two quartic coupling constants not 
zero, see Fig. |l]b, in contrast to the tricritical O(N) point where the 
quartic term disappears and a sextic term takes over |L9 . 



c. N < N c : F_ and F + move into the complex parameter space, see fig. |T]c 
and d. The absence of stable fixed points is interpreted as a signature 
of a first order transition. 

The difficulty is to find a reliable value for N c . The dependence on d has 
been calculated to first order in e = 4 — d ITBL ITSfl as 



N c = 21.8 - 23.4 e . (2) 

In repeating the analysis for the triangular antiferromagnet and choosing 
e = 1 for d = 3 Kawamura has argued that the XY and Heisenberg cases 
should be above N c , motivated by Monte Carlo calculations which supported 
a second order transition of a new universality class || f|, |5], ||, |7|, |8| with 
critical exponents for the Heisenberg case v = 0.59, (3 = 0.29 (see table [l]) 
different from the standard ones. 

Of course the results linear or higher in e are at best asymptotic [f20| 



Antonenko, Sokolov and Varnashev extended the analysis to e 2 M. For an 
analysis directly in three dimensions see |K|. For the critical number of 
components after resummation 

N c = 3.91(1) (3) 
is obtained by the direct method [H| and N c = 3.39 with the e expansion 



to second order 0. In Fig. |2] the results of the RG studies are depicted. 
There is a line which separates a region of first order transition for high 
dimensions and low iV (Fig. [l]c) from a region of second order transition for 
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low dimensions and high iV (Fig. |T]a). The transition on this line is special 
(Fig. |I]b) but is not the standard tricritical. 

Similar considerations have been applied to the normal to superconduct- 
ing type II transition. Starting point is the Ginzburg-Landau Hamiltonian 
with also two coupling constants, a quartic coupling and the coupling to the 
magnetic field. Expanding to first order in e = 4 — d Halperin et al. |21j and 
Chen et al. [22| arrived at results almost identical to the ones for frustrated 



magnets. There are also four fixed points with one stable one if the number 
of order parameters N exceeds N c = 365.9 for a fictitious superconductor in 
4 dimensions. The unstable 0(2iV) fixed point is replaced by an XY fixed 
point. For N < N c only the unstable Gaussian and this unstable XY fixed 
point are present, and the question whether the physical relevant N = 2 case 
is of second order or not can also not be answered, although the value of N c 
in two loop order has been calculated 
If we put e = 1 we obtain N c (d = 3) 
transition should be of second order 



N c = 365.90 - 640.76 e + 0(e 2 ) . 
—275 which is less than two and the 
The nematic-to-smectic-A transition 



in liquid crystals is described by a model similar to the Ginzburg-Landau 
Hamiltonian transition [24 , [25f and again the problem of second or first order 
transition arises [21, p6[ . 

The superfluid phase transitions of He 3 and the phase transitions of he- 
limagnets were first studied with the field theoretical approach ([I]) by Love, 
More, Jones and Bailin [16|, [17| and by Garel and Pfeuty |18| respectively. 
In principle all the considerations should also be applicable to these systems, 
however for He 3 the critical region is not really accessible. Returning to the 
frustrated antiferromagnets and accepting the value N c = 3.91 for d = 3 of 
eq. the XY and Heisenberg systems and also the helimagnets should have 
a first order transition in three dimensions. 



2.2 Almost second order transitions and complex fixed 
points 

The discrepancy between the results of RG studies which indicate a first 
order transition and the results of Monte Carlo study which favor a second 
order transition (see table [|) can be removed by using complex fixed points 
in the RG analysis. Indeed, as N c is very close to N = 3 the fixed points 
F + and F_ will not be far off from the real space of parameters u* and v* at 
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N c as depicted in fig. |T]d. Zumbach [|TJj] has studied in detail the influence 
of a complex fixed point on the RG analysis. The basin of attraction of such 
a complex fixed point will mimic a second order transition with modified 
scaling relations 



2(3 = u(d-2 + r]-cp) 
7 = v (2 - r] + c 7 ) . 



(4) 
(5) 



The constants cp and c 7 are corrections of the scaling relations absent if the 
fixed point is real. Using these relations to determine r\ without corrections 



we get a negative value for r\. Since rj must be positive jO] or at least zero we 



have an estimate for cp a and an indication that the transition is not a real 
second order one. We have used this criterion before for x-y case [[y|. The 
essential point is the smallness of rj, approximately r) ~ 0.03 for magnets, so 
that the correction terms cp }J are clearly visible. 

A negative rj is unphysical since unitarity would be violated in the corre- 
sponding quantum field theory |11], [TJ, p^j . For spin glasses this restriction 
no longer holds and indeed a negative rj can be obtained |28||. Also for su- 



perconductors a negative rj have been found |22| influenced by the choice 
of gauge [p9|1 . The critical exponents must be gauge-independent and in a 
recent study JHJ it has been shown that r\ is positive for superconductors. 

Using the local potential approximation for the RG the critical exponent 
v can be calculated for real and complex fixed points (more precisely for a 
minimum in the flow) following Zumbach 
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In this approximation rj is 
zero and the result will have errors of a few percents, for example 8% for 
the Heisenberg ferromagnet |]3T|. For the fixed point F + of the frustrated 
systems Zumbach obtains v ~ 0.63 for iV = 3. This result is compatible with 
v = 0.59(1) obtained by MC calculations for Heisenberg spins on a stacked 
triangular lattice (see table |]). 

By the same approximation Zumbach estimated the interval where the 
systems should be under the influence of a complex fixed point as 



2.58 < iV < 4.7 



(6) 



Outside of these limits for iV < 2.58 the transition is of first order and for 
iV > 4.7 the transition belongs to the new universality class. The two values 
iV c = 4.7 for the local potential approximation and N c = 3.91(1) found by 
Antonenko and Sokolov |10[ by field theoretical means are not so different. 
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The minimum limit 2.58 is in agreement with our simulations for the Stiefel 
model. The transition is indeed of first order for N = 2 [13] and it will be 
shown here that the N = 3 case has a pseudo second order behavior. 



3 Models and simulations 

In this section we introduce briefly different models that we use in this work. 
A more complete presentation can be found in |0|. 



3.1 The STAR model 

For the stacked triangular antiferromagnet (STA) we take the simplest Hamil- 
tonian with one exchange interaction constant J > (antiferromagnetic). 

H = H^i S *- S i ' ( 7 ) 

(y) 

where Sj is a three component spin vector of unit length, and the sum is over 
all neighbor spin pairs of the lattice. There are six nearest neighbor spins 
in the plane and two in adjacent planes. In the ground state the spins are 
in a planar arrangement with the three spins at the corners of each triangle 
forming a 120° structure (see Fig. 

For one triangle the spin vectors obey the equation 

Si + S 2 + S 3 = . (8) 

where the indices refer to the three corners. To get the STAR model one im- 
poses this equation as a restriction for the spin directions valid at all temper- 
atures. In this theory the local fluctuations violating this constraint become 
modes with a gap. Thus they do not contribute to the critical behavior and 
can be neglected. However, one can do this only for selected triangles, for 
instance the shaded ones in fig. |3], otherwise the spin configuration would be 
completely frozen. 

The lattice is partitioned into interacting triangles which do not have 
common corners, so that each spin belongs only to one triangle. For the 
Monte Carlo simulation one spin direction (two degrees of freedom) is cho- 
sen and then the direction of the second spin selected in the cone of 120° 
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around the first one. The third is determined by the constraint (|8]). All the 
spins of the supertriangles (the shaded ones) have, only in the ground state, 
the same orientation. At finite temperatures there is a perfect order within 
a supertriangle, but fluctuations of the orientation between these triangles 
occur. 

The Monte Carlo updating for the state of the supertriangle is done the 
following way. First two orthogonal unit vectors are chosen at random in 
three dimensions. One needs three Euler angles to do this, the first 9q must 
be chosen with probability sin(6 l o) d9o and the other two 9%^ with probability 
6*1,2 so that 

/ sin(6 l o) cos( 

u = sin(0 o ) sin(0i) | (9) 
cos 





/ — SI) ! ( II i j (■•( )s( {I-, ) — COS (0o ) cos l 

cos(0 o ) sin(0i) sm(0 2 ) | (10) 
\ sm(0 o ) sin(0 2 ) 

Then the three spin vectors on each supertriangle are determined by 

51 = u 

5 2 = cos(120°)u + sin(120°)v (11) 

5 3 = cos(240°) u + sin(240°) v . 

The interaction energy between the spins of this supertriangle with the 
spins of the neighboring ones is calculated in the usual way. We follow the 
standard Metropolis algorithm to update one supertriangle after the other. 
One million MC-steps for equilibration are carried out and up to six million 
steps were used for the largest sizes to obtain reliable averages. Long enough 
simulations are necessary because the critical slowing down is strong. 

The number of spins is given by N = L x L x L z , where L x L gives 
the number of spins in one plane and L z = 2L/3 the number of planes. 
Simulations have been done for L = 18, 21, 24, 30, 36, 42, where L must be 
a multiple of 3 in order to use periodic boundary conditions and to avoid 
frustration effects in the planes. 

The order parameter M used in the calculations is 

M = ^EI^I (12) 

i=i 
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where the magnetization Mj is defined on one of the three sublattices. This 
definition generalizes the one used for two collinear sublattices. 



3.2 The Stiefel model 

The Stiefel model V^, that is a "Zweibein" (or dihedral model) in three di- 
mensional space is a further abstraction of the constrained three spin system 
discussed in the previous section. The three spins at the corner of a triangle 
can be taken as a planar or degenerate "Dreibein" where the third leg is a 
linear combination of the other two and can be left out. The energy is given 
by 

H = JJ2[ eiW • e 1 {j) + e 2 « • e 2 (j) } (13) 

ij 

where the mutual orthogonal three component unit vectors ei(i) and e 2 (i) at 
lattice site i interact with the next pair of vectors at the neighboring sites j. 
The interaction constant is here negative to favor alignment of the vectors at 
different sites. Also a cubic lattice instead of a hexagonal one can be taken 
(for further details see fl"3H). 



We use a single Monte-Carlo cluster algorithm j32] . A cluster of connected 
spins is constructed and updated in the standard way. The first site of the 
cluster o is chosen at random together with a random reflection r. Then all 
neighbors j are visited and added to the cluster with probability 

P = 1 - exp(min{0, E x + E 2 }) (14) 

where 

£ 1>2 = 2/3(r-e li2 ( ))(r-e 1)2 (j)) (15) 

and j3 = J/T until the process stops. Finally all spins of the cluster are 
"reflected" with respect to the plane _L to r, that is 

eJ2" = eg-2(eg.r)r. (16) 

It has been demonstrated that this cluster method is very efficient in reducing 



the critical slowing down for the O(N) ferromagnet p3 |. We thought it 
should be useful in our case since the order is also ferromagnetic, contrary 
to the original STA model where the spins have a 120° structure so that the 
cluster algorithm described does not work. However, we did not observe a 
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decreasing of the critical slowing down. Indeed on general grounds Sokal et 



al. |34| have argued that Wolff's cluster algorithm cannot reduce the critical 
slowing down unless the manifold for the order parameter is a sphere as for 
O(N) ferromagnet. 

In each simulation 5 million measurements were made after enough single 
cluster updatings of 1 million steps for equilibration. Cubic systems with 
linear dimension L = 15, 20, 25, 30, 35, 40 were simulated. In comparison 
with the STAR model L should be multiplied by V3 since one site of the 
Stiefel model represents three spins. The equivalent sizes would be 25 to 70. 

The order parameter M for this model is 

M=^g|A/,| (17) 

where Mj is the total magnetization given by the sum of the vectors e, over 
all sites and N = L 3 is the total number. 

3.3 The right handed trihedral model 

If one adds a third vector given by the vector product 

e 3 (i) = ei(z) x e 2 (i) . (18) 
and takes an energy formally identical to the previous one 

H = J£>i(i)-ei0') + e 2 (z)-e 2 (j) + e 3 (z) • e 3 (j) ] (19) 

(y) 

has one constructed a really different model? Since introducing this new vec- 
tor e 3 no degree of freedom has been added and therefore this system should 
belong to the same universality class as the original V 3j 2 Stiefel model. It is 
not the V 3j3 Stiefel model since only right handed systems are constructed. 
An additional Ising symmetry is absent |1| contrary to the x-y or planar 



case where the "chirality" as an additional Ising variable is present. 



In the work of Azaria et al. |35|, [3(| it has been shown that the interaction 
between the third vector, which is absent in the original STA model, will be 
automatically generated by the renormalization group in d = 2 + e expansion 
and, at the fixed point, should be equal to the first two. 
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This system has been studied by Loison and Diep |37j before. Here we 
will confirm that the transition is strongly first order. The Monte Carlo 
procedure is similar to the one used for the STAR model. In addition to the 
two orthogonal vectors u and v of eq. |9|-|T0| a third vector given by the vector 



product w = u x v is constructed. The interaction energy between the spins 
of this trihedral with the spins of the neighboring trihedral is calculated in 
the usual way and we follow the standard Metropolis algorithm to update 
one trihedral after the other. In each simulation 10 000 Monte Carlo steps 
were made for equilibration and for taking the averages. Cubic system of 
linear dimension up to L = 30 were simulated. The order parameter M is 
similar to the previous case 

1 3 

M= — y\MA (20) 
i=i 

but one has to take the average of three Mj instead of two. 



4 Numerical method 

4.1 Finite size scaling for second order transitions 

We use the histogram MC technique developed by Ferrenberg and Swendsen 



[ |38|j and divide the energy range into 30,000 intervals, for more details see 



H| . The errors are determined with the help of the Jackknife procedure [|3G 

For each temperature T we calculate the internal energy per site E = 
(E)/N and the specific heat C = ((E 2 ) - (E) 2 ) / \N T 2 ) , where (. . .) indicate 
the average and iV is the number of sites. Similarly we determine the aver- 
ages of the order parameter or staggered magnetization M = (M) and the 
corresponding susceptibility x — N ((M 2 ) — (M) 2 )/T. The quantities needed 
besides M for the finite size analysis are defined below 

X2 = j(M 2 ) (21) 
N 3 

X4 = ^((M 4 > -3(M 2 ) 2 ) (22) 

U - l-J^L (23) 
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(24) 



where \2 is the magnetic susceptibility and Xi the fourth derivative of the 
free energy with respect to the magnetic field in the high temperature region 
where the order parameter is zero. The cumulant V\ is used to obtain the 
critical exponent u, and the fourth order cumulant U to determine the critical 
temperature. 

According to the FSS theory [40, for a second order transition the 
various quantities just defined should scale for a sufficiently large system at 
the critical temperature T c as 

(25) 
(26) 
(27) 
(28) 

where (3 = 1/T and g are constants not dependent on size L. We will not 
use x to determine j/u but xi (^) and % 4 fl2"6"|) using 

(29) 



X2 


= 9x2 


L ~//u 


X4 


= 9 X a 




Vi 


= 9Vi 


L l/u 


M 


= 9m 





74/1/ = d + 27/z/ 



since the errors are smaller. 

To find the critical temperature T c we record the variation of U with T 
for various system sizes and then locate T c as the intersection of these curves 
[42]j , since the ratio of U for two different lattice sizes L and V = bL should 
be 1 at T c , that is 



Ur 



1 . 



(30) 



T=T C 



Due to the presence of residual corrections to finite size scaling, actually one 
has to extrapolate the results taking the limit (lnfe) -1 — > (Fig. [5|, |7| and 



4.2 First order transitions 

A first order transition has a different scaling behavior ^5 
a. The histogram P(E) should have a double peak. 
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b. Magnetization and energy should show a hysteresis. 



c. The minimum of the fourth order energy cumulant W = 1—(E 4 ) / (3(E 2 ) 2 ) 
varies as 

W = W* + bL~ d and W* ^ § . (31) 

A double peak in P(E) means that at least two states with different energies 
coexist at the same temperature. As a consequence the fourth order energy 
cumulant cannot be |. This fact was employed for the smaller sizes in a pre- 



liminary study ||37|| . Hysteresis effects have to be expected in the simulation 
of larger systems where the two peaks are well separated since the transi- 
tion time from one state to the other grows exponentially with the system 
size. For frustrated antiferromagnets we are studying these criteria are not 
very helpful. Only for the trihedral model the hysteresis is clearly visible in 
Fig. [l(| indicating a strong first order transition. 



5 Results 

5.1 STAR model 

Using (^) we first determine T c . U is plotted for different sizes from L — 18 
to L = 24 as a function of temperature in Fig. [|. From the intersections we 
extrapolate T c as 

T c = 1.43122(12) , (32) 

see Fig. The estimate for the universal quantity U* at the critical temper- 
ature is 

U* = 0.6269(10). (33) 

With the value of T c we determine the critical exponents by log-log fits. 
We obtain v from V\ (Fig. ^|), 7/1/ and 74/V from % 2 and Xi (Fig- |; n °t 
shown), and j3/u from M (Fig. [J): 



V 


= 0.504(10) 


(34) 


l/v 


= 2.131(13) 


(35) 


74/ v 


= 7.252(29) 


(36) 




= 0.439(8) . 


(37) 
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The uncertainty of T c is included in the estimation of the errors. The value of 
7/V = 2.126(14) found from 74/^ using (|29"D is compatible with the one found 
directly. Combining these results we obtain (3 = 0.221(9), 7 = 1.074(29) and 
from the scaling relation 

7/^ = 2-77 (38) 

we obtain 77 = —0.131(13). The results are summarized in table f|. 

A negative value of 77 is impossible, for a second order phase transition 
it should always be positive [[□]]. A slightly negative value was already dis- 
cernible in some of the results of the unconstrained STA-model, see table 
|l|. In making use of our previous analysis of the x-y case |H| we conclude 
that the critical behavior should be described by the renormalization flow 
of a fixed point in the complex parameter space. The scaling relation for 77 
above is then modified to give a positive value for this critical exponent, see 
previous section. Moreover this large negative value cannot come from only 
the presence of a complex fixed point. Thus we interpret the negative value 
as a crossover from the basin of attraction of the complex fixed point to the 



true first order transition, see Fig. [XT 



In a previous study of the STAR model, Dobry and Diep |46j obtained for 
the exponent v = 0.440(20) which disagrees with ours 0.504(10). However 
this seems to be a misprint since the inverse value is also given and v = 
1/2.08 = 0.48 is within the statistical errors. However the other exponents 
seem wrong. We repeated these calculations for larger clusters and with 
better statistics in order to have evidence that this system shows a weak first 
order transition, since unambiguously 77 is really negative 



5.2 Stiefel's V^ 2 model 

Following the same procedure as for the STAR model we determine first T c . 
The cumulant U plotted as a function of temperature for different system 
sizes from L = 15 to L = 40 is shown in Fig.f7| The extrapolation, shown 
in Fig. U, for L = 20 and L = 25 agrees rather well; the little difference for 
L — 15 indicates that this lattice size is not yet sufficient for a strictly linear 
extrapolation. We obtain 

T c = 1.5312(1) (39) 
and for the universal quantity U at T c 

U* = 0.6326(12). (40) 
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In plotting the logarithm of Vi, X2, Xa arid M as function of the logarithm of 
temperature difference (T — T c )/T c , shown in Fig. |9], one finds from the slope 

v = 0.507(8) (41) 

j/v = 2.240(10) (42) 

lA / v = 7.480(22) (43) 

P/u = 0.381(6). (44) 

From 74 /V using (|2~9"D we get j/u = 2.240(11) which is the same value found 
directly from \2- Further we obtain for (3 = 0.193(4) and 7 = 1.136(23) and 



with the scaling relation (pq) an even more negative value 77 = —0.240(10). 
See for a listing table |, where also the specific heat exponent a = 2 — dv has 
been added. The errors given include the uncertainty of the estimate of T c . 
The results of Kunz and Zumbach Jl5| are also given there and are com- 



patible with ours. Their a is actually determined by fitting specific heat data 
in the high temperature region and v with the hyperscaling relation. They 
noticed also that 77 is negative, apparently their value —0.10(5) is too small, 
compared to —0.24(1) we found. We cannot share their conclusion that the 
negative 77 is due to a strong finite size effect not properly taken care of. 
We rather take it as an indication that an analysis as a genuine second or- 
der transition leads to inconsistencies and interpret the negative value as a 
crossover from the basin of attraction of the complex fixed point to the true 
first order transition, see Fig. [□]. 

5.3 Direct trihedral model 

It is known from |37j] that the energy distribution P(E) near the critical 
temperature has a double peak structure. A rough estimate of the correlation 
length £ by h of the smallest size where the two peaks are well separated 
by a region of zero probability gives £0 = 6. Here we determine the energy 
cumulant at the transition. Extrapolating W to a system of infinite size we 
obtain 

W* = 0.623(1) (45) 

which of course deviates from the value | for a second order transition. It is 
not difficult to see hysteresis effects. In Fig. [H] hysteresis for the energy E 
is clearly visible. 
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Therefore the behavior of a system of trihedrals differs from the one of 
dihedrals. In weakening the interaction between the third components of 
the trihedrals one could of course recover the apparent second order nature 
of the dihedral system. The first order character of the transition of pure 
dihedrals is impossible to show directly. This system remains always under 
the influence of the virtual or complex fixed point for the sizes of systems 
one can simulate, while the presence of the third vector allows the system 
of trihedrals to stay outside of the neighborhood of this fixed point and the 
"true" first order behavior is seen directly. 

6 Discussion 
6.1 Summary 

The STAR and the Stiefel model have slightly different critical expo- 
nents. Looking at table ^, most noticeable are the differences for f3 and 7: 
Pstar = 0.221(9) and (3 sue = 0.193(4), outside the estimated statistical and 
systematic errors. The difference to the values of (3 ~ 0.286 for the orig- 
inal stacked triangular antiferromagnet in table [I] is even larger. Also the 
non-collinear antiferromagnet on the bet lattice has a different /3-value (see 
table 1). Using the scaling relation rj = 2 — d + 2(3/v the exponent 77 is al- 
ways negative. Moreover, a member of the same universality class, the right 
handed trihedral model, has a strong first order transition. This does not 
correspond to the standard behavior of a universality class one expects from 
the seemingly clear evidence for a second order transition as visible in Fig. |9] 
for the STAR and the Stiefel model. 

In taking the point of view that the simulations show an extremely weak 
first order transition it is very natural to use a field theoretical description 
with the usual renormalization group scheme supplemented by complex fixed 
points, see Fig. [I]. Such complex fixed points exist if the number of compo- 
nents N is lower than a critical one N c , and according to field theoretical 
analysis the critical dimension is N c ~ 4. As for a real fixed point the 
renormalization flow will be attracted and the system imitates a second or- 
der behavior, however it can escape if the system size L is larger than the 
largest correlation length £0 possible and then the crossover to first order 
region is reached. This crossover phenomenon is known to occur in the two 
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dimensional Potts model with q = 5 components H7L H8 



In the analysis based on complex fixed points Zumbach [14| showed that 
the scaling relations should be modified by correction terms eqs. (Q-^)- Leav- 
ing out these terms the scaling relations for rj produces negative values, for 
XY spins of the STA model it is slightly negative (rj ~ —0.06 fll3| ) and for 
the models studied here the effect is two or three times larger, see table [|. 
Another possibility of obtaining negative values for rj is to visualize the spin 
system in the crossover region between a second to first order transition. Us- 
ing the same scaling relations as before the exponent 77 will tend to the value 
f] = — 1 of a "weak first order transition" (table |3|). Our MC simulations of 
the STAR and the Stiefel model are in such a crossover region otherwise the 
very negative values for rj « —0.2 cannot be understood. The same applies 
for MC simulations of helimagnetics on bet lattices ||. 

In contrast to the XY case with N = 2 the frustrated Heisenberg case 
is closer to the critical dimension. That is, the STAR and the Stiefel model 
appear at first sight as second order transitions whereas in the XY case they 
show really first order behavior [TJ] . The stacked triangular antiferromagnet 
STA follows the same trend with a small negative 7] in the XY case and an 
only slightly negative one in the Heisenberg case table |l|. The smaller sphere 
of influence or basin of attraction of the complex fixed point in the XY case 



compared to the Heisenberg case is depicted in Fig. O and Fig. 112 



Another interpretation following Kawamura |49| takes the diagram in 
Fig. |l]a as basis, where an attractive F + and an unstable fixed point F_ 
are neighbors. The initial points in the RG flow for the STAR and the 
Stiefel model are placed outside the domain of attraction of the fixed point 
F + , that is, to the left of the line (G,F_) in Fig. |I|a while the STA model 
should be to the right of this line. Therefore the STAR and the dihedral 
model have a first order transition and the STA model a standard second 
order transition. The difficulty is the critical dimension N c > 3 according 
to Antonenko and Sokolov [[IIJ and actually Fig [Tjc should be the starting 
point. Beyond all question is that all the MC simulations give a negative 
rj impossible to explain with a real fixed point JTT|. One observes also the 



experimental crossover behavior from second order to first order in systems 
belonging to the XY class like the holmium, dysprosium or CsCuCl 3 (see 
[ Q| ) which is in disfavor of this interpretation. 
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6.2 Experiments 



In real systems because of the omnipresence of planar or axial anisotropies 
the behavior of frustrated Heisenberg antiferromagnets will be difficult to 
observe. Moreover such systems are usually quasi one or two dimensional, 
that is a succession of crossovers will occur from 2d to 3d and from Heisenberg 
to Ising or XY behavior. To observe the first order behavior the temperature 
must be very close to the critical temperature in order to notice that the 
correlation length is limited in the basin of attraction of the complex fixed 
point F. For XY spin systems, in spite of this limitation the crossover has 
been observed in certain materials like Ho fl50l , Dy [5^| , CuCoCl 3 [^3 
(see also 
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Before the first order region is reached in Heisenberg systems the crossover 
to Ising or XY behavior prevents a first order transition of Heisenberg type. 
In Fig. [13] a typical scheme is drawn for the crossovers and transitions of a 
system with Ising anisotropy, for more details see ffl and for other crossovers 



3- 

Nevertheless the second order transition can be studied. A list of results is 
in Table |: for VC1 2 ]55[, for VBr 2 @, for Cu(HCOO) 2 2CO(ND 2 ) 2 2D 2 |7 
and for Fe[S 2 CN(C 2 H 5 ) 2 ] 2 Cl ||58| . For the last two examples the observed ex- 
ponents might be influenced by the crossover between 2d to 3d Ising behavior. 
The experimental results agree quite well with the MC simulation in table |l[ 



6.3 Renormalization group studies in 2+e 

For the STA model studied in d = 2 + e expansion using the corresponding 
nonlinear a-model Azaria et al. f55| , |3"E|| obtained for all iV a stable fixed point 
F 2+e at a small distance e from d = 2. In comparing results in lowest order 
of 1/N for the 2 + e to the 4 — e expansion they show that the fixed points 
found in the different approaches are equivalent for large N. 

Specially for iV = 3 a clear evidence by symmetry arguments can be 
given that the transition near d = 2 should be of the standard ferromagnetic 
0(4) class. If the transition is of second order also for d = 3 it must therefore 
belong to the 0(4) class with the exponent 704 = 1.54 (table |3|) very different 
from the value of the STA model ^sta — 1-18 found by the MC method 
(table |T]). This is not compatible with the Monte Carlo simulations and 
experimental studies (see tables). 
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We think that we can rule out the possibility of a ferromagnetic 0(4) 
class for N = 3 if the relevant operators are the same for large and low N. 
The 1/N expansion for 7 (similar arguments hold for all exponents) for both 
the ferromagnetic non frustrated (NF) |60| case and frustrated case (F) [0 
are: 



JNF 



If 



d-2 
2 

d-2 



1 



6 

N 
9 



where 



siii(7r(d-2)/2)r(tZ-l) 



(46) 
(47) 

(48) 



2 7rr(rf/2) 2 

is a positive quantity. Thus for iV large 7a?f > If- The reason is the 
different breakdown of symmetry from 0(N) to 0(N — P) with P = 1 in 
the ferromagnetic case and P = 2 in the frustrated case. When N decreases 
the ratio P/N increases also the difference between 7^ and 7^ will increase 
and one should have ^nf > If for any N and in particular for N = 3. 

However it is possible that irrelevant operators in the 4 — e and 1/N 
expansions become relevant near dimension two and low N. In this case 
for iV = 3 the transition should be of first order for dimension three and a 
little less than three. Below this limit the transition could be of the "chiral" 
universality class and in approaching two dimensions it should become a 
4 transition. This hypothesis has found some support f5T|. We note that 



at exactly two dimensions, systems are known to be driven by the 0(4) 



symmetry, at least at low temperature f62 ]. 
Other RG studies near two dimensions 



have found more than 



a single transition governed by the 0(4) fixed point with an intermediate 
"nematic" phase between the disordered high temperature region and the 
ordered 120° structure at low temperatures. This double transition should 
also occur for the right handed trihedral model studied here when the inter- 
action between the third components are larger than between the two others. 
The ordering of the third vector occurs before the ordering of the other two 
vectors and both transitions are of second order, the first of Heisenberg type 
and the second of XY type. The intermediate "nematic" phase can also be 
interpreted as a ferromagnetic phase and the two transitions characterized 
as paramagnetic-collinear and collinear-noncollinear. 
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Another possibility discussed in the literature is that the transition is 
influenced by the presence of topological defects which are not visible in the 
continuum and perturbation theory as in RG (2 + e, 4 — e and directly in 
d = 3) f62|| . In this work we have assumed that the effects of topological 
defects are irrelevant for the critical behavior Hl5 



7 Conclusion 

We have tried to discuss the phase transition of frustrated Heisenberg spin 
systems in general terms. The starting point is Monte Carlo simulations of 
systems for which the condition of local rigidity is imposed. From the finite 
size analysis we suggest that the transition is of first order. Our result is 
that one can rely on the renormalization group studies and the true behavior 
is indeed first order. However one should use the concept of complex fixed 
points to describe the almost second order behavior of frustrated Heisenberg 
spin systems. In the strict sense there is no "chiral universality class" but 
in practical terms it exists. We have given reasons why in the usual Monte 
Carlo simulations and in experiments this first order transition is difficult or 
even impossible to observe. 
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Table 1: Critical exponents by Monte Carlo for 0(3). STA = stacked 
triangular antiferromagnet, bet = body-centered-tetragonal, calculated by 
!) 1 j y = 2 - rj , 2 ) dv = 2 - a , 3 > 2 (5jv = d - 2 + rj . 
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Table 2: Experimental values of critical exponents for compound VCI2, 
VBr 2 , 1 Cu(HCOO) 2 2CO(ND 2 ) 2 2D 2 and 2 Fe[S 2 CN(C 2 H 5 ) 2 ] 2 Cl. 
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Figure 1: a), b), c) Hamilto- 
nian flow induced by renormalization 
group transformations. The arrows 
indicate the direction of flow under 
iterations, d) For N < N c the fixed 
points F + and F_ become complex. 
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Figure 3: Ground state con- 

figuration for the STA and the 
STAR model. The supertriangles are 
shaded. The dihedral are drawn at 
the center of each supertriangle. 
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Figure 2: The critical curve N c (d) 
separating the first order region at 
small N and large d from the second 
order region for large N and low d. 
Curve (1) is the result of first order 
in e = d — 4, and curve (2) in sec- 
ond order. The STA at d = N = 3 
is marked by a cross. In contrast to 
curve (1) curve (2) lets the physical28 
relevant systems iV = 2, 3 stay inside 
the first order region. 



Figure 4: The interaction between 
the "spins" or dihedrals % and j of 
Stiefel's model V 3t2 '■ where ei(i) in- 
teracts only with e\(j) and e 2 (i) with 
e 2 (j)- 



0.68 



0.660 
U 



0.640 



0.620 



0.600 




1.428 



1.430 



1.432 



Figure 5: Binder's parameter U for 
the STAR model as function of the 
temperature for different sizes L (in 
the left part of the figure upwards 
from L—18 to L=42). The arrow in- 
dicates the critical temperature T c . 
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Figure 6: Estimated T c for the STAR 
model plotted vs inverse logarithm of 
the scale factor b = L'/L . For clarity, 
only the results for L=18, 21,24 are 
shown. The estimated temperature 
is T c =1.43122(12). 
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Figure 7: Binder's parameter U for 



Stiefel's model V-. 



3,2 



as function of 



temperature for different sizes L (in 
the left part of the figure upwards 
from L = 15 to L = 40). The arrow 
indicates the critical temperature T c . 
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Figure 8: Estimated T c for the Stiefel 
model : plotted vs inverse loga- 
rithm of the scale factor b = V / L . 
For clarity, only the results for L = 
15, 20, 25 are shown. The estimated 
temperature is T c —1. 5312(1). 
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Figure 9: Values of X2, V\ and M as 
function of L in ln-ln scale at T c for 
the STAR and the Stiefel model : 
to get the slopes 7/V, 1/V and 
Error bars are smaller than symbols. 
To obtain a clear figure values for M 
are shifted by +3 units. 
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Figure 10: Hysteresis in the internal 
energy per spin E versus T for the 
right handed trihedral system. Lines 
are guides to the eye. The arrows in- 
dicate if the MC simulation is cooling 
(circle) or heating (square) the sys- 
tem. 
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Figure 11: Hypothe- 
sis on the Hamiltonian flow induced 
by renormalization-group transfor- 
mations for N < N c (see Fig. |T] c) 
showing the smaller influence of the 
complex fixed point for the XY spins 
(N= 2) compared to the Heisenberg 
spins (iV=3). The arrows give the di- 
rection of flow and the black circles 
indicate the positions accessible for 
the different systems in simulations 
or experiments. The hatched region 
is under the influence of the complex 
fixed point and the outer circle repre- 
sents the crossover to the first order 



region. 
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Figure 13: Crossover phenomena in 
experimental systems (not to scale). 
There is a crossover to Ising behav- 
ior before reaching the first order re- 
gion. F + indicates that the behavior 
is caused by the complex fixed point. 
For real systems it is possible to have 
other crossover behaviors due to lat- 
tice distortions and Dzyaloshinsky- 
Moriya interactions. 
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